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ABSTRACT 


A STUDY OF THE SHOCK WAVE STRUCTURE IN 
UNDEREXPANDED JETS EXHAUSTING INTO STILL AIR 

Kenneth Michael Larsen 
September 1991 

The main objective of this project was to predict the expansion wave/oblique shock 
wave structure in an underexpanded jet expanding from a convergent nozzle. The shock 
structure was predicted by combining the calculated curvature of the free pressure 
boundary with principles and governing equations relating to oblique shock wave and 
expansion wave interaction. The procedure was then continued until the shock pattern 
repeated itself. A mathematical model was then formulated and written in Fortran to 
calculate the oblique shock/expansion wave structure within the jet. 

In order to study shock waves in expanding jets, Schlieren photography, a form of 
flow visualization, was employed. Thirty-six Schlieren photographs of jets from both a 
straight and 15° nozzle were taken. An iterative procedure was developed to calculate the 
shock structure within the jet and predict the non-dimensional values of Prandtl primary 
wavelength (w/r n ), distance to Mach Disc (L<j) and Mach Disc radius (r^). These values 
were then compared to measurements taken from Schlieren photographs and experimental 
results published by NASA in Refs(3) and (13). 

The results of the method discussed above agreed closely to measurements from 
Schlieren photographs and data from Refs(3) and (13). This method provides excellent 
results for pressure ratios below that at which a Mach Disc first forms. Calculated values 
of non-dimensional distance to the Mach Disc (Lj) agreed closely to values measured from 
Schlieren photographs and published data. The calculated values of non-dimensional Mach 
Disc radius (rj), however, deviated from published data by as much as 25% at certain 
pressure ratios. 
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CHAPTER 1 
Review of Literature 


Jet Structure 

When a gas exhausts from a nozzle into still air it will undergo an expansion or 
compression depending on the nozzle exit static pressure to ambient pressure ratio (pj/p 0 ). 
If this ratio is greater than one, the nozzle is underexpanded, and if the ratio is less than 
one, it is overexpanded Ref(l). In the underexpanded case the jet exhausting into still air 
will expand in order to equalize with ambient pressure. Expansion or compression at the 
nozzle exit can occur whether the nozzle contour is convergent (known as sonic) or 
convergent-divergent (known as supersonic). Expansion or compression of an exhaust jet 
can be accomplished for given nozzle inlet conditions by variation of the nozzle back 
pressure. 

In this project, a convergent nozzle mounted on a six component thrust stand was 
supplied by compressed air and exhausted into ambient air. Schlieren photographs were 
taken at various inlet pressures for the cases of a straight and 15° convergent nozzle. In 
each case photographed, the nozzle was underexpanded. In addition, the nozzle total 
pressure to ambient pressure ratios exceeded the critical pressure ratio in all cases 
photographed and consequently the flow was choked (M = 1.0) Ref(2). 

The critical pressure ratio (P n /pj) for air behaving as an ideal gas with a constant 
specific heat ratio of 1.4 is 0.52828 Ref(2). Therefore, the minimum nozzle total pressure 
required to choke the flow for an ambient air pressure of 14.7 psia is 27.8 psia. It is of 
interest to note that the shock diamond structure is first observed in Schlieren Photograph 
No. 22 and is provided in Figure 1. The corresponding nozzle total pressure in this 


photograph is 31.7 psia. Since shock waves cannot form in subsonic flow, the nozzle must 
be choked before shocks can appear in the jet. 

As the jet expands from the nozzle, a constant pressure boundary forms Ref(3). 
This boundary indicates the region in which the static pressure of the jet equals ambient 
pressure. The slope of this boundary is a maximum at the nozzle exit (where expansion is 
greatest) and gradually diminishes until the jet stops expanding and compression begins 
Ref(4). The result is a parabolic-shaped boundary with a point of inflection where 
compression begins. 











Figure 2 illustrates the theoretical structure of the jet expanding from an 
underexpanded convergent nozzle into still air Ref(l), (31) and (43). 



Figure 2. Jet Structure of an Axisymmetric Jet Exhausting Into Still Air 

Since the flow is choked, the Mach number at the exit of the nozzle (Region 1) is 1.0. At 
this point the flow is assumed to expand to atmospheric pressure isentropically via Prandtl- 
Meyer expansion Ref(3). As the flow expands from Region 1 to 2, it accelerates to 
supersonic speed and its static pressure decreases to atmospheric pressure. In addition, 
static temperature and density decrease as the flow expands. 

The expansion from Region 1 to 2 is modeled as a discrete expansion wave. In 
actuality, there is an expansion fan consisting of an infinite number of Mach waves Ref(3). 
This discrete expansion wave is the average of the initial and final Mach waves. The 
calculation of the angle of the discrete expansion wave with respect to the flow axis is 
shown in Appendix A. From symmetry, the angle of the discrete expansion waves from 
the upper and lower edges of the nozzle are equal. Between Regions 1 and 4, the discrete 
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expansion waves intersect and reflect at the centerline of the flow. Expansion waves of 
equal strength reflect in the same manner as reflection from a solid surface. Therefore, as 
the flow passes from Region 2 to 4, it encounters another expansion wave. As a result, its 
static pressure decreases to a value less than atmospheric pressure. 

The expansion waves then reflect from the free pressure boundary. From the 
principles of expansion waves and oblique shock waves (Appendix A), expansion waves 
reflect from free pressure boundaries as oblique shock waves. Therefore, from Region 4 
to 5 in Figure 2, the flow passes through a reflected oblique shock wave. Since the flow 
passes through an oblique shock wave its static pressure increases to atmospheric pressure 
in Region 5 in order to equalize with the free pressure boundary. The flow Mach number 
decreases and static pressure, temperature and density increase as the flow passes through 
an oblique shock wave Ref(2). The flow structure on the upper and lower sections of the 
jet in Figure 2 are assumed to be symmetric. 

The oblique shock waves then intersect at the centerline between Regions 4 and 7. 
Oblique shock waves of equal strength reflect in the same manner as a reflection from a 
solid surface. Therefore, as the flow passes from Region 5 to 7, it encounters another 
oblique shock wave. The static pressure in Region 7 is therefore greater than atmospheric 
pressure. Finally, the reflected oblique shock waves meet the free pressure boundary and 
reflect as expansion waves Ref(2). 

The oblique shock waves in Region 7 reflect from the constant pressure jet 
boundary as expansion waves. These expansion waves would then intersect and reflect at 
the centerline in an identical fashion to the flow from Region 1 to 4. The flow pattern 
beyond Region 7 therefore repeats itself. 
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The point at which the jet shock pattern repeats itself is known at the Prandtl 
primary wavelength (w/r n ) and the repeating shock diamond structure was first detected by 
L. Prandtl using Schlieren photographs Ref(5). The shock structure discussed and 
illustrated in Figure 2 occurs at pressure ratios below that at which a Mach Disc first 
appears. The Mach Disc is a normal shock wave that appears between Regions 4 and 7. 
The Prandtl primary wavelength is an important aspect of the study of axisymmetric free 
jets. This point is readily visible in all Schlieren photographs. The length (w/r,j) is a non- 
dimensional distance and is easily measured from Schlieren photographs. Results of 
measurements of Prandtl primary wavelength measured from photographs taken in this 
project are compared with results from other references in the results section of the report. 


Constant Pressure let Boundary 

The first step in the calculation of the jet structure was to determine the contour of 
the constant pressure jet boundary at various nozzle pressures. The structure of the 
expansion/oblique shock waves within the jet was then calculated using the governing 
equations for oblique shocks, expansion waves, and normal shock waves in conjunction 
with the jet boundary contour. The jet boundary contour was provided in the form of a 3rd 
order polynomial. 

The constant pressure jet boundary contour was solved by using the method of 
characteristics. Analytical work was conducted by D. C. Pack in 1947 Ref(6) and 
compared with Schlieren photographs in 1959 by Love, Grisby, Lee and Woodling Ref(3). 
L. Prandtl first used Schlieren photographs to study jets exhausting form a convergent 
nozzle or orifice and was discussed in Ref(5). Love and Grisby used the method of 
characteristics to determine the contour of the constant pressure jet boundary. In this 
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report, an alternative method of predicting the constant pressure jet boundary contour 
proposed by Adamson and Nicholls Ref(4) was used. 

The first step in the solution of the jet boundary was to determine the initial flow 
angle of the jet at the nozzle exit. If the air was assumed to undergo a Prandtl-Meyer 
expansion at the nozzle exit then the initial flow angle was given as follows Ref(2). 

a = u 2 = -(b)^ tan-1 [ (M 2 2 - 1.0)/b ]' A + tan- 1 (M 2 2 - 1.0) 1/4 Eqn(l) 

7 + 1 

where: b = 

7-1 

The derivation of Equation 1 is provided in Appendix A. As discussed earlier, this was the 
maximum slope of the constant pressure jet boundary. 

The flow angle or tangent of the jet boundary curve at any point is simply the 
difference between the initial flow angle (a) and the Prandtl-Meyer angle (5) at that 
location. The slope of jet boundary at any point can therefore be put in terms of the 
vertical and horizontal coordinators r and x as follows. 


tan (a - 6) = dr/dx 


Separating variables, integrating and putting r and x into non-dimensional form r/r n and 
x/r n where r n is the nozzle radius. 


x/r n = 


fr/r n 


1 


d(r/r n ) / tan(a - S) 


Eqn(2) 
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The integration of Equation 2 can be solved numerically using either the Newton-Cotes, 
Simpson's 1/3 or Simpson's 3/8 methods Ref(ll). 

Adamson's paper Ref(4) suggested that the integration method should give good 
results for relatively low pressure ratios compared with the method of characteristics. The 
results of jet boundary calculations using the Adamson-Nicholls method at pressure ratios 
of 2, 10, and 20 are compared with results from Love, Grisby Ref(3) in Figure 3. This 
method provides an accurate prediction of the jet boundary contour at low pressure ratios 
and departs from the method of characteristics solution as the pressure ratio increases. 
More importantly, the calculated boundary contour closely approximates the actual 
boundary measured from Schlieren photographs. Since the pressure ratios encountered in 
this project range from 1.1 to 3.2, the Adamson-Nicholls method of Ref(3) was sufficient. 
A convenient feature of the use of non-dimensional variables was that different experiments 
could be compared, independent of the nozzle dimensions used in each experiment. 


Oblique/Expansion Wave Structure Within a Jet 

The next area of interest in the study of jets exhausting into still air was the 
formation and interaction of expansion and shock waves in the flow. Previous work in this 
area used the method of characteristics to predict the jet structure Ref(3). However, the 
method of characteristics was limited to those jet structures without Mach Discs. In this 
report, the prediction of the jet structure was simplified when reflection and intersection of 
discrete expansion waves, normal and oblique shock waves in quasi-one-dimensional flow 
were considered. The principles of expansion and oblique shock wave interaction are 
presented in Appendix A. 
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The general structure of the axisym metric jet has been described. With this general 
structure known, in addition to the nozzle inlet conditions, the governing equations of 
expansion waves and oblique shock waves were used to determine the jet structure as a 
function of inlet pressure. That is, with nozzle inlet conditions known, the governing 
equations of oblique shock waves and Prandtl-Meyer expansion were used to determine the 
expansion wave angles 0), the oblique shock wave angles (e) and the flow turning angles 
( 5 ) throughout the flow structure. With these values known in addition to the curvature of 
the jet boundary, the jet structure was defined. Consequently, the Prandtl primary 
wavelength (w/r^ was calculated for a given nozzle inlet pressure. The exact calculational 
procedure for predicting the jet structure was provided in Appendix B. 

Love, Grisby Ref(3) developed empirical relations (based on measurements from 
Schlieren photographs) for predicting the Prandtl primary wavelength for the case of no 
Mach Disc present and for the presence of a Mach Disc. The two equations are as follows. 

No Mach Disc present (Pj/p 0 < 2.0): 

w/r n = 3. 1[ Mj 2 (pj/ Po ) - 1 ]' A 

Mach Disc present (pj/p 0 > 2.0): 


w/r n = 3.04(pj/p o ) 0 - 437 + 3.1 [ (2 M^ - \) l/l - 1 ] - (0.55) ^ 
+ 0.5 { (1/1.55) [ ( Pj/po - 2 ) ] >/4 - 1 } 


Where: 


w/r n = Non-dimensional Prandtl primary wavelength 
Mj = Mach number at nozzle exit (Mj = 1.0 in this case) 
/9j = (Mj - 1.0)' 1/4 = 0.0 (for Mj = 1.0) 
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For the case of a convergent nozzle (Mj = 1.0) these equations simplify as follows. 


No Mach Disc present: 


w/r n = 3.1 [ (pj/po) - 1 ]' A 


Eqn(3) 


Mach Disc present: 


w/r n = 3.04(pj/p o ) 0 437 - 0.5 


Eqn(4) 


Formation of the Mach Disc 

The next step in the study of jet structures was to model the Mach Disc 
phenomenon. The case illustrated in Figure 2 is known as a regular oblique shock wave 
reflection. In this case, the reflected oblique shock between Regions 5 and 7 is strong 
enough to straighten out the flow to become parallel to the original flow axis Ref(l). A 
Mach Disc forms when a regular reflected oblique shock wave from Region 5 to 7 cannot 
straighten out the flow. At pressure ratios above this point, a Mach Disc becomes visible 
in Schlieren photographs. This phenomenon of irregular oblique shock wave reflection is 
known as Mach reflection Ref(7). Figure 4 illustrates the change in shock diamond 
structure from Region 4 to 7 when a Mach Disc forms. 



Figure 4. Mach Disc Structure 

The Mach Disc forms when the initial oblique shock wave formed between Regions 
4 and 5 generates a reflected oblique shock with a corresponding deflection angle greater 
than the maximum deflection angle (& maT ). That is, when the deflection angle 
corresponding to the Mach number in Region 5 (behind the incident oblique shock) is 
insufficient to straighten out the flow, a Mach Disc forms Ref(6). The maximum 
deflection angle associated with an oblique shock wave is given as follows and was derived 
in Appendix A. 

r - 1 pi 

& = Where: r = Eqn(5) 

2(r) 2 p\ 

Additional information desired with regards to the Mach Disc is the non- 
dimensional radius Od) and non-dimensional distance to the Mach Disc (Lj). The Mach 
Disc is modeled as a normal shock wave. In actuality, it has a slight curvature Ref(7). 
A Mach Disc is a normal shock wave that connects four oblique shock waves. In Figure 4, 
Regions 4, 5 and 6 coincide with those in Figure 2. 


In the case of a Mach Disc another area of interest appears, the region behind the 
normal shock. From Figure 4, the Mach numbers in Regions 7, 8, and 9 are different. 
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The only flow property in common between Regions 7, 8, and 9 is the static pressure. 
Since there is a difference in properties between these regions a discontinuity appears in the 
form of a slip-line Ref(8). The slip-line is visible in some photographs in Ref(7), 
however, it does not appear in any photographs taken in this project. Higher resolution 
photography should visualize the slip-lines. 

The next step in the study of Mach Disc formation was to determine the size and 
location of the Mach Disc. That is, the non-dimensional distance to the Mach Disc (Lj) 
and the non-dimensional radius Od). Analysis involving the method of characteristics was 
unable to accurately predict the geometry of the Mach Disc. Recall the Love, Grisby 
empirical relationship for the Prandtl primary wavelength with a Mach Disc present in the 
jet (Equation 4). 

w/r n = 3.04(pj/ Po ) 0 437 - 0.5 

If w/r n is known by actual measurement or by use of Equation 4, then the structure 
of the Mach Disc can be determined using data calculated with the Jet Boundary and Jet 
Structure programs. A calculational procedure for predicting the Mach Disc distance and 
radius as a function of the nozzle pressure was developed and is presented in Chapter 3 and 
a sample calculation provided in Appendix G. 
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CHAPTER 2 

Method of Experimentation 
Principles of Schlieren Photography 

An important aspect of this project was the flow visualization of the air discharging 
from the straight and 15° convergent nozzles. The technique employed was Schlieren 
photography. The method was first developed by Foucault (1859) and Toepler (1864) 
Ref(10). L. Prandtl first used Schlieren photography to study the shock wave phenomenon 
in axisymmetric jets. The basic principle behind Schlieren photography is that certain 
types of discontinuities or property variations in a transparent media refract light. The 
result is an image or trace of the discontinuity onto a photographic plate in the form of a 
shadow. In this case air was exhausted from a convergent nozzle at high pressure into still 
air. The air undergoes changes in static density, pressure and temperature as it flows out 
of the nozzle and expands. In addition, expansion waves, normal, and oblique shock 
waves form in the flow. As light passes through the discontinuity or shock waves, the 
light is refracted resulting in an image of the shock which can then be photographed. 

The primary components of the Schlieren apparatus used in this project and 
illustrated in Figure 5 are the light source, light slit, two Schlieren mirrors, knife edge and 
photographic plate (35 mm camera in this case). This variation of Schlieren system is 
known as the Barnes and Bellinger, or two mirror system Ref(10). The light source used 
was a Halogen lamp placed at a distance corresponding to the focal length of the first 
mirror. Each Schlieren mirror has a 48 inch focal length and a 4 inch diameter. The 
maximum field of view of a photograph is therefore limited to 4 inches. There is a slit 
located between the light source and first Schlieren mirror. The light slit concentrates the 
light onto the first Schlieren mirror. The Schlieren mirrors are positioned so that the flow 
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field is located between them. A knife edge is placed at the focal length of the second 
mirror. The purpose of the knife edge is to reduce the intensity of the light exposing the 
photographic surface and create a dark ground condition, Ref(10). The light source, slit, 
and knife edge are placed parallel to one another. 

Schlieren photographs were taken with and without the use of the knife edge. 
Photographs made without the knife edge gave greater detail of the flow. The photographs 
presented in Appendix D were made at a wide range of nozzle total pressures. 
Photographs were also taken using two different techniques of camera arrangement. In the 
first method, light passes through the knife edge which shines directly onto a 35 mm 
camera with the lens removed. This method required the use of low ASA film (125) and 
fast shutter speeds (1/500 to 1/1000 sec). The absence of a knife edge under these 
circumstances results in severe overexposure of the film regardless of the ASA used. In 
the second method, the beam of light was reflected from a flat reflective surface and a 
photograph was taken of the reflection. This method allowed photographs to be taken with 
or without a knife edge. Unlike the first method, this method required the use of 
extremely high ASA film (3200) and low shutter speeds (1/4 to 1/8 sec). 

The best results were achieved using the second method of camera arrangement. 
The Schlieren system is conceptually simple and was set up easily. However, the difficulty 
of Schlieren photography occurred when experimenting with shutter speeds and film type 
to achieve the best photographic results. In addition, complete darkness in the laboratory 
is essential to ensure proper film exposure. 
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Schlieren Mirror #2 
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Thrust Stand Data Acquisition 

The primary component of the thrust stand assembly are the convergent nozzle 
(straight and 15°), air inlet plenum, thrust stand, isolation mounts or bellows, load cells 
for force measurement and pressure transducers. One pressure transducer is located in the 
inlet plenum and measures static pressure. The second pressure transducer is located just 
before the nozzle inlet and measures total pressure. There were 6 load cells measuring 
forces in the x, y, and z directions. 

Data from the pressure transducers was sampled by the Fluke Data Acquisition 
System where the date was displayed and printed. A printout of pressure and force 
measurements was made for each Schlieren photograph presented and an example is 
provided in Appendix J. Each printout contains the nozzle inlet total pressure and load cell 
readings from each of the six component load cells. Additional experimental data that 
needed to be determined was the system total temperature (storage tank temperature) and 
barometric reading. 
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CHAPTER 3 
Analytical Analysis 


Recall that in the Constant Pressure Jet Boundary section of Chapter 1 that the 
integral representing the jet boundary defined by Equation 2 was as follows: 

f r/r n 

x/r n = d(r/r n ) / tan (a - 6) Eqn(2) 

J l 

Where 

5 = Difference in Prandtl-Meyer angles at any location along the jet boundary. 

S = u 2 - vi = -(b )' 1/4 {tan 1 [(M 2 2 - 1.0)/b] l/4 - tan^fM, 2 - 1.0)/b] ,/4 } 

+ [tan 1 (M 2 2 - 1.0)' 1/6 - tan- 1 (M, 2 - 1.0)'^] 

q = Initial jet turning angle defined by Equation 1. This value remained constant 
throughout the calculation and was a function only of the nozzle pressure (Pn). 

a = i/ 2 = -(b) 1/4 tan 1 [ (M 2 2 - 1.0)/b ]' A + tan 1 (M 2 2 - 1.0) 1/4 Eqn(l) 

b= . 

7-1 

The integral in Equation 2 can be numerically calculated. The Simpson's 3/8 rule of 
integration was defined as follows Ref(l 1). 

x = (3/8) H (F 0 + 3F a + 3F b + 2F C + 3F d + 3F e + 2F f + 3F g -I- 3F h + Fj) 

H = Amount that the r coordinate is incremented (variable). 


x = Distance along the flow axis. 


The numerical procedure is outlined as follows. 
Evaluate function of F 0 : 
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The value of F 0 was the value of the function within the integral sign in Equation 2. 

F 0 = 1/tan (a - So) 

The value S 0 was the difference in Prandtl-Meyer angles at the start of the integration. In 
this case it was the nozzle exit so (6 0 = 0). F c then becomes. 

F 0 = 1/tan (a) 

Evaluate function values F a through F^ 


1. Select a value for H and increment r by H. 

r a = r n + H 

2. Evaluate the critical area ratio Ref(2) at point a. 

(A/A*) a = (r a / rn ) 2 (A/A*)e 

A* = Nozzle exit area (0. 1625 in 2 ) 
r n = Nozzle radius (0.2275 in) 

(A/A*)e = Critical area ratio associated with the Mach number at the 
end of the Prandtl-Meyer expansion at the nozzle exit. 

Solve for (A/A*) E 

<a/a-) e 

Where Me is the Mach number after the initial Prandtl-Meyer expansion from the 
nozzle and was defined by Equation 10a in Appendix A. Substitute r a , r n and 
(A/A*)e into the above equation and solve for (A/A*) a . 


r + 1 


( 1 + 7-1 M e 2 ) i 

A J 


7 + L 
2 ( 7 - 1 ) 
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3. 


Solve for Mach number at location a (MJ implicitly from the following equation. 


(A/A*) 

a 


_L_ 


M a 2 


i ± i 
2 ( 7 - 1 ) 

r 2 ( 1 + i - 1_ M a 2 ) i 

7+1 2 J 


4. Solve for « a by substituting M a and Me into the expression for 6 defined above. 

5. Solve for the function value at location a (Fa)- 

F a = 1/tan (a - Sj 

6. Increment r a by H 
r b = r a + H 

7. Return to Step 2 and solve for (A/A*) b 
Where: (A/A*) b = 0b/r n ) 2 (A/A*) E 

8. Calculate F b . 

F b = 1/tan (a - 6 b ) 

9. Repeat procedure until Fj is calculated. 

10. Substitute values of F 0 through Fj into the Simpson's 3/8 Rule equation and solve 
for x. Make values of r and x non-dimensional by dividing by the nozzle radius r n 
and store values. This completes the integration for one data point. 

11. Set F 0 = Fj and return to Step 1 to begin integration again. 

12. Repeat the procedure until the desired number of jet boundary data points are 
generated. 


The above procedure for calculating the jet boundary was programmed into the Jet 
Boundary Fortran program presented in Appendix F. The program generated values of 
r/r n and x/r n along the jet boundary in accordance with the Adamson-Nicholls integration 
technique of Ref(4). The only data that the user must input into the program was the 
nozzle total pressure (P n ). Values of r/r n and x/r n from the Jet Boundary program were 
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then put into a curve fitting program in order to calculate a 3rd order polynomial 
representing the jet boundary. This process was completed for pressures corresponding to 
the Schlieren photographs. The resulting jet boundary contours are presented in 
Appendix E. 

Calculation of .Tet Structure (No Mach Disc Present) 

The general jet structure illustrated in Figure 2 and discussed in the Jet Structure 
section of Chapter 1 represented the case of no Mach Disc present in the flow. For given 
nozzle inlet conditions, the discrete expansion wave angles 0*), oblique shock wave 
angles (c) and flow turning angles (5) throughout the flow can be calculated using the 
governing equations of oblique shock waves and Prandtl-Meyer expansion. The jet 
boundary contour was calculated using the Adamson-Nicholls integration method and a 3rd 
order polynomial representing the curve was generated. With the expansion/oblique shock 
wave angles known in addition to the boundary contour, the geometry of the jet structure 
can be determined. Figure 2 is repeated at this time for convenience. 



Figure 2. Jet Structure of an Axisymmetric Jet Exhausting Into Still Air 
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The calculation^ procedure for determining the jet structure and ultimately the Prandtl 
primary wavelength (w/rjj) for the case of no Mach Disc present in the flow is presented 
below. 

1. Calculate the Mach number in Region 2 with nozzle total pressure (PJ known. 
Substitute into Equation 10a of Appendix A. 

2. Solve for the initial flow turning angle (a) using Equation 1. 

3. Solve for the angle of the discrete expansion wave (average Mach wave angle /J ave )- 
This was the angle the discrete expansion wave at the nozzle exit makes with 
respect to the flow axis. 

4. Solve for the discrete expansion wave angle (with respect to the flow axis) from 
Region 2 to 4. This was the expansion wave reflected at the flow centerline. 

5. Introduce the jet boundary contour 3rd order polynomial and solve for non- 
dimensional distance (xi + X2). 

6. Using the governing equations of oblique shock waves, the geometry of the shock 
structure and the jet boundary contour, solve for non-dimensional distance 

(X3 + X4). 

7. Determine the Prandtl primary wavelength w/r n by adding the non-dimensional 
distances xi through X4. 

w/r n = xi + X2 + X3 + X4 


This calculational procedure was carried out in the sample calculation provided in 
Appendix B. The Jet Structure program was written to calculate the Prandtl primary 
wavelength (w/r n ). The program requires the user to input nozzle total pressure (P n ), total 
temperature (T n ) and the four coefficients associated with the jet boundary 3rd order 
polynomial. The program is presented in Appendix G. 
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Calculation of .Tet Structure (Mach Disc Present) 


Figure 4 is repeated at this time for convenience. 



Figure 4. Mach Disc Structure 


The procedure for calculating the non-dimensional distance to the Mach Disc (Lj) and 

Mach Disc radius (r^) is as follows. 

1. Solve for flow properties in Regions 1 through 6 in the same manner as the case of 
no Mach Disc present. 

2. Solve for the non-dimensional distances xj, X2 and X3 in the same manner as in the 
case of no Mach Disc present. 

Define: x/r n = Xi + X2 + X3 

3. Use the governing equations of oblique shock waves, solve for the oblique shock 
wave angle £45 and flow turning angle 54s for the initial reflected oblique shock 
from Region 4 to 5. 

4. Using the governing equations of normal shock waves, solve for the flow properties 
in Region 9. 

5. Set the static pressures in Regions 7, 8, and 9 equal. 
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6. Solve for the flow properties in Regions 7, 8, and 9 and the oblique shock wave 
angle £57 and the flow turning angle 657 (reflected oblique shock). 

7. Verify that « 57 is less than the 5,^ corresponding to the Mach number in Region 5 
using Equation 5. 

8. Determine the value of the Prandtl primary wavelength w/r n from either Equation 4 
or measurement from a Schlieren photograph. Use w/r n , the value of x/r n 
determined in Step 2, the jet contour and the geometry of the Mach Disc structure 
in Figure 4 to calculate values of r d and L d . 

Where r d = Mach Disc radius 

L d = Distance to Mach Disc 

9. Compare calculated values of L d and r d with those determined from curve fitting 
NASA graphical data from Ref(10). These equations were determined as follows. 


r d = -0.948 + 0.673 (pj/p 0 ) - 0.086 (p/p 0 ) 2 + 0.0043 (pj/p 0 ) 3 


L d = 0.055 + 1.396 (pj/p 0 ) - 0.198 (pj/po) 2 + 0.013 (pj/p 0 ) 3 

A sample calculation of the above procedure was conducted for the case of 
P n = 89.0 psia and is presented in Appendix C. The Jet Structure Fortran program 
was written to carry out this calculation and is presented in Appendix G. 
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CHAPTER 4 
Experimental Results 


Schlieren Photographs 

Thirty-six Schlieren photographs were taken of the straight and 15° convergent 
nozzles. These photographs are presented in Appendix D. Figure 6 (Photo No. 29) is an 
example of a Schlieren photograph of the straight nozzle at a pressure ratio at which there 
is no Mach Disc present in the flow. Figure 7 (Photo No. 30) is also a Schlieren 
photograph also of a straight nozzle. However, in this photograph a Mach Disc is present 
in the flow. 

Values of Prandtl primary wavelength (w/r^ measured from Schlieren photographs 
taken in this project for the straight and 15° nozzles are provided in Table 1. In each case 
the pressure ratio is below that at which a Mach Disc first forms in the flow. These values 
are compared with values predicted by the NASA empirical relation Equation 3 from 
Ref(3). 
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Figure 6. Schlieren Photo No. 29 
P n = 75.6 psia 
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Figure 7. Schlieren Photo No* 30 
P n = 88*6 psia 






Table 1 

Comparison of Prandtl Primary Wavelength Data (No Mach Disc) 
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Photo No. 

P„ (psia) 

Pj/Po 

w/r n (photo) 

w/r n Eqn(3) 

15° Nozzle 

13 

35.4 

1.27 

1.62 

1.63 

14 

42.1 

1.51 

2.10 

2.21 

15 

52.0 

1.87 

3.00 

2.88 

16 

58.4 

2.10 

3.29 

3.25 

17 

65.9 

2.37 

3.50 

3.62 

18 

68.2 

2.45 

3.88 

3.73 

Straight Nozzle 

23 

35.5 

1.20 

1.62 

1.60 

24 

42.0 

1.51 

2.14 

2.21 

25 

51.7 

1.86 

2.68 

2.87 

26 

58.3 

2.10 

3.33 

3.24 

27 

66.0 

2.37 

3.71 

3.63 

28 

67.8 

2.44 

3.71 

3.71 

29 

75.7 

2.72 

4.00 

4.07 


In Table 1, measurements from Schlieren photographs agree very closely with those 
determined by NASA in Ref(3). 


Recall that in the Constant Pressure Jet Boundary section of Chapter 1 , the first step 
in the solution of the jet boundary was to determine the initial flow angle (a) of the jet. In 
Appendix A, the procedure for calculating this angle was presented. The following table 
compares the calculated values of initial flow angle with those predicted in Ref(3). 





31 

Table 2 

Initial Jet Flow Angle 


Photo No. P n (psia) Pj/p 0 <*° a Ref(3) 



In the sample calculation provided in Appendix B the pressure ratio was (pj/p 0 = 
2.0). The Prandtl primary wavelength (w/r n ) was calculated as 3.16. Substituting pj/p Q 
into Equation 3, the value of NASA's predicted value of w/r n was 3. 1 . The two values are 
within 2%. Values of calculated Prandtl primary wavelength (w/r n ) predicted by NASA 
Equation 3 and values measured from Schlieren photographs in this project are compared 
with those predicted using the Jet Boundary and Jet Structure programs. Results are 
summarized in the following table. 


Table 3 

Comparison of Prandtl Primary Wavelength Data With Jet Strucutre Program 

Results (No Mach Disc) 
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The Jet Structure program predicts Prandtl primary wavelength (w/r n ) values quite close to 
those predicted by NASA Equation 3 and values obtained from Schlieren photograph 
measurements. 

Calculation of Mach Disc Formation Pressure 

In the sample calculation in Appendix B, pj/p 0 was 2.0, corresponding to a nozzle 
total pressure of 55.65 psia. The value for (£)max in Region 5 to 7 was calculated as 14.5° 
which was greater than the (5) required to straighten out the flow (12.5°). Therefore, no 
Mach Disc was formed at this pressure ratio. 

A piece of information that is desired is the total pressure or static pressure ratio at 
which is first exceeded in Region 5 to 7 (reflected oblique shock) and consequently 

the pressure at which a Mach Disc should form. This can be calculated using the Jet 
Boundary and Jet Structure programs. The following calculational procedure is presented. 

1. Increment the nozzle total pressure above 55.65 psia and run the Jet Boundary 
program to calculate the constant pressure boundary. 

2. Insert jet boundary data into a curve fitting program to generate a 3rd order 
polynomial representing the jet boundary. 

3. Input above data into the Jet Structure program to see if (djjnax is exceeded. If not, 
increment the total pressure and return to Step 1. 

The above procedure was completed and the nozzle total pressure at which the 
Mach Disc should first appear was calculated as 72.5 psia. This prediction was 
substantiated by experimental work conducted by Addy in Ref(12) and data provided by 
NASA in Ref(13). Addy's paper predicts that the Mach Disc will first form at a nozzle 
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total pressure ratio (P n /Po) of 5.0 or (P n = 73.5 psia). In Schlieren Photo No. 28 
(P n = 68.4 psia) there is no indication of a Mach Disc, in Photo No. 29 (P n = 75.6 psia) 
there is a very small Mach Disc and Photo No. 30 (P n = 88.6 psia) has a distinct Mach 
Disc. Additional Schlieren photographs in this pressure range, with use of a higher quality 
light source, should further substantiate this prediction. The close agreement between 
Addy Ref(12), Schlieren photos 28 through 30 and predictions from the Jet Boundary and 
Jet Structure Fortran programs further validate the procedure for shock structure prediction 
presented in this report. 

.Tet Structure Calculation (Mach Disc Present) 

A sample calculation of the Mach Disc structure was provided in Appendix C. This 
procedure was also programmed into the Jet Structure program. The value of r d 
represented the non-dimensional radius of the Mach Disc. From Schlieren Photograph No. 
30 this value was (0.45). The Mach Disc radius was plotted as a function of nozzle static 
pressure ratio in Ref(13). The value is 0.42 for pj/p 0 = 3.2 (P n = 89.0 psia). The value 
of L d from the sample calculation in Appendix C was 3.007. This value was also 
measured from Schlieren Photo No. 30 and was 3.0. This also agrees with the NASA data 
in Ref(13). 

Values for Mach Disc distance (L^) and radius (r^) plotted as functions of nozzle 
static pressure ratio in Ref(13) were provided in Appendix H. Data from these graphs 
were curve-fit, resulting in the following equations. 

Non-dimensional distance to Mach Disc (L^: 


L d = 0.055 + 1.396 (pj/ Po ) - 0.198 (pj/ Po ) 2 + 0.013 (pj/ Po )3 


Eqn(6) 
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Non-dimensional Mach Disc radius Od): 

r d = -0.948 + 0.673 (pj/p 0 ) - 0.086 (pj/po) 2 + 0.0043(pj/p o ) 3 Eqn(7) 

Curve fitting enables this empirical data to be programmed into the Jet Structure program. 
The following Table compares L d and r d data as calculated by the Jet Structure program, 
Equation 6 and Equation 7, and measurement from Schlieren photographs. 


Table 4 

Values of (L d ) and (r d ) as a Function of Pressure Ratio 


Photo No. Pn 

Pj/Po 

Ld(prog) Ld(Eqn6) 

%Dev 

Ld(photo) 

rd(prog) 

rd(Eqn 7) 

% Dev 



72.5 

2.61 

2.81 

2.58 

8.2 






— 

29 

75.7 

2.72 

2.85 

2.65 

6.9 

2.8 

— 

— 

— 

— 

80.0 

2.87 

2.90 

2.75 

5.4 

— 

0.394 

0.378 

4.1 

— 

85.0 

3.05 

2.96 

2.85 

5.0 

— 

0.403 

0.428 

5.8 

30 

89.0 

3.20 

3.01 

2.93 

2.7 

2.95 

0.402 

0.466 

13.5 

31 

95.0 

3.41 

3.08 

3.04 

1.3 

3.0 

0.403 

0.514 

21.8 

— 

100.0 

3.59 

3.14 

3.13 

0.4 

— 

0.416 

0.555 

25.0 

— 

110.0 

3.95 

3.28 

3.29 

0.5 

— 

0.451 

0.593 

22.1 

— 

115.0 

4.11 

3.35 

3.37 

0.6 

— 

0.513 

0.663 

23.7 

.... 

120.0 

4.29 

3.43 

3.45 

0.6 


0.531 

0.695 

23.5 


From Table 4, the values of Mach Disc distance (L^ from NASA predictions, Schlieren 
photographs and the Jet Structure program, agree closely over a wide range of pressure 
ratios. Values of Mach Disc radius (r^), however, deviate as much as 25% from NASA 
experimental data. 


Conservation of Mass Across the Mach Disc Structure 

The flow properties across the oblique/normal shock wave structure in Regions 4, 
5, and 9 were calculated using the Jet Structure program. The dimensions of the Mach 
Disc structure were calculated in the sample calculation in Appendix C for P n = 89.0 psia. 
Recall the Mach Disc structure. 
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Figure 4a. Mach Disc Structure 

Where: A4J = surface area of the oblique shock wave 

Ao = vertical component of A45 
A9 = surface area of the normal shock wave 

The conservation of mass across an oblique shock wave is expressed as follows Ref(2). 

P4 V 04 = ps Vnj Eqn(8) 

V n = velocity normal to the oblique shock wave 

Normal velocity is related to free stream velocity as follows. 


V04 = V 4 sin (£45) 

Where (£45) is the oblique shock wave angle. Expressing V n in terms of Mach number. 

Vn4 = M4 a4 sin (£45) = M4 (7 R t4) ,/4 sin (£45) 

The normal velocity on the other side of the shock wave (Vos) is related to the free stream 
velocity in Region 5 as follows. 

Vos = M5 a 5 sin (£45 - S45) = M 5 (7 R t 5 ) ,/4 sin (£45 - £45) 

Equation 8 can be rewritten as follows. 

p 4 M4 (7 R U)‘ A sin (£45) = ps M 5 (7 R t 5 ) ,/5 sin (£45 - £45) 

Writing the continuity equation across the Mach Disc. 


rru = ms + 1TI9 
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Making substitutions. 

p* M 4 (7 R U)'^ A 4 = p$ Ms (7 R ts ) ,/4 sin (£45 - $ 45 ) + P9 M 9 (7 R t9) A A 9 Eqn(9) 

From Figure 6a, A45 is realted to A 0 as follows. 

A 45 = A J sin (e 4 s) 

Also from Figure 6a, A 0 can be calcualted from the known non-dimensional values r ra and 
r d- 

A 0 = it ( r ra (0.2275) ) 2 = * ( r d (0.2275) ) 2 

In the above expression r ra and r d are made dimensional by multiplying by the nozzle 
radius (r^. For the case of P n = 89.0 psia. 

r ra = 1.355 
r d = 0.402 

Substituting values and solving for A 0 , A 4 J and A 9 (ms = 43.2°). 

A 0 = 0.2725 in 2 
A45 = 0.3981 in 2 
A 9 = 0.026 in 2 

The values of static density, temperature, pressure and Mach number for Regions 4, 5, and 
9 were calculated using the Jet Structure program and are as follows. 

pt = 0.0488 lbm/ft 3 P9 = 0.173 

M 4 = 2.69 M 9 = 0.496 

t4 = 214.5°R t9 = 500.3 

ps = 0.1184 
M 5 = 1.682 
t 5 = 335.5 

After substituting values into Equation 9 the left and right hand sides of the equation 
balance. 

0.1956 = 0.1956 


And mass is conserved across the Mach Disc. 
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Discussion 

Schlieren photograph measurements of Prandtl primary wavelength (w/r n ), Mach 
Disc radius (r^ and Mach Disc distance (Lj) consistendy agreed with measurements 
conducted by NASA and presented in Refs(3) and (13). This can be observed by 
examination of Tables (1) and (4). In addition, litde deviation was observed between 
measurements of Prandtl primary wavelength from photographs of the straight and 15° 
nozzles at similar pressure ratios. 

The oblique shock/expansion wave structure within the jet was calculated by 
combining the calculated jet boundary contour in the form of a 3rd order polynomial with 
principles and governing equations of oblique shock/expansion wave interaction and the 
basic geometry of the jet as illustrated in Figure 2. Sample calculations were provided in 
Appendix B and C and the Jet Boundary and Jet Structure Fortran programs were written 
to calculate the jet structures with and without Mach Discs. Results of calculations are 
presented in Tables 3 and 4. 

Figure 3 indicated that the use of the Adamson-Nicholls integral method for 
calculation of the jet boundary contour resulted in contours close to those determined by 
the use of the Method of Characteristics. The deviation between the two methods at a 
pressure ratio of pj/p 0 = 2.0 was approximately 4%. Since the pressure ratios encountered 
in this project were in the range of 1.1 to 3.2, the Adamson-Nicholls method was a 
sufficient model for jet boundary contours. 

The Jet Boundary and Jet Structure programs were combined in order to predict the 
pressure ratio at which the Mach Disc would first appear in the flow and was calculated as 
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72.5 psia. This prediction is substantiated by Addy Ref(12), NASA Ref(13) and by the 
examination of Schlieren photographs taken in this project. 

Table 3 summarized the Jet Structure program results for predicted Prandtl primary 
wavelength (w/r^ for the case of no Mach Disc present in the jet. The Table also 
compared these values with those predicted by NASA in Ref(3). The results agree quite 

closely. Table 4 compared predicted values of Mach Disc distance (Lj) and Mach Disc 

radius Od) using the Jet Structure Program and those determined experimentally by NASA 
Ref(13). Calculated values of Mach Disc distance (Lj) are consistently in close agreement 
with NASA empirical data. Analytical prediction of Mach Disc radius, however, deviated 
as much as 25% from empirical data. Conservation of Mass across the Mach Disc 

structure was investigated in order to validate flow properties calculated by the Jet 

Structure program. 


Conclusions 


Schlieren photography provided an indispensable source of information for 
comparison with analytical models presented in this report for the prediction of the 
axisymmetric jet expansion/oblique shock wave structure. However, Schlieren 
photographs were unable to visualize certain aspects of the flow structure including slip- 
lines associated with the Mach Disc structure. In addition, there was poor definition of 
small Mach Discs, particularly in the nozzle total pressure range of 75 psia to 80 psia. 
These deficiencies can be avoided in future Schlieren photography conducted on the thrust 
stand with use of a higher quality Schlieren system light source. 

The Adamson-Nicholls integral method of computing the contour of the constant 
pressure jet boundary agreed closely with the Method of Characteristic solution at low 
pressure ratios. If larger pressure ratios are encountered in subsequent work with the 
thrust stand, then the accuracy of the calculated jet boundary contour could degrade, 
resulting in a less precise jet structure approximation. A possible remedy for this problem 
would be to incorporate a Method of Characteristics solution into the Jet Boundary 
program. This would allow the two methods to be compared at any supply pressure. 

The jet structure calculation method discussed and illustrated by example in 
Appendix B produced excellent results at pressure ratios at which no Mach Disc formed in 
the flow. Although the expansion from the nozzle exit was modeled as Prandtl-Meyer 
expansion (isentropic flow/perfect gas) and the resulting expansion fans were approximated 
as discrete expansion waves, excellent results were obtained for the Prandtl primary 
wavelength (w/r n ) estimates. 
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Calculation of the jet structure with a Mach Disc present required several 
simplifying assumptions. Recall that in the sample calculation in Appendix C, Equations 
30 through 34 provided 5 equations with 7 unknowns. The simplifying assumptions made 
were that the oblique shock wave angles £45 and t$i were equal and that the Prandtl 
primary wavelength (w/rj was known. Resulting estimates of Mach Disc distance (L<j) 
were accurate, however, calculated values of Mach Disc radius were not as accurate. The 
method overall functions well as a means of approximating the Mach Disc structure. No 
previous analytical work of this kind was found in the review of literature. 
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PrandtI-Mever Expansion/Initial Jet Boundary Angle (a) 

When the air initially expands to atmospheric pressure at the nozzle exit in the 
underexpanded converging nozzle, it is assumed to encounter a Prandtl-Meyer expansion. 
In Prandtl-Meyer expansion air behaves as a perfect gas expanding isentropically from an 
original Mach number of 1.0 or greater to a region of lower static pressure and 
consequently higher Mach number. This is illustrated in Figure 8. 



Figure 8. Expansion of a Sonic Flow 

In Figure 8, (^0 is the Mach line angle associated with Region 1 (90° for M = 1.0), fa) 
is the Mach line angle in Region 2 and (£) is the flow turning angle. The governing 
equation for Prandtl-Meyer expansion is as follows, Ref(4). 

6 = V 2 - 1 /, = -(b) 1/4 { tan- 1 [ (M 2 2 - 1.0)/b ]* - tan-' [ (M, 2 - 1.0)/b ]* } 

+ [ tan-' (M 2 2 - 1.0)* - tan-' (M, 2 - 1.0)* ] 

In this equation (vi) is the Prandtl-Meyer angle in Region 1 (0° for M = 1.0), fa) is the 
Prandtl-Meyer angle in Region 2. In the case of Mi = 1.0, the equation simplifies as 


follows. 
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a = vi = -(b) 1/4 tan-i [(m 2 2- 1.0)/b] l/6 + tan'(M 2 2 - 1.0 )' ' h 

Where: 


b = 


i + 1 


7-1 


Eqn(l) 


In Equation 1, vi is replaced by a to avoid confusion in subsequent calculations. This 
equation gives the initial turning angle of the flow at the nozzle exit. As discussed earlier, 
this is the maximum slope of the constant pressure jet boundary. There are two unknowns 
in Equation 1, M 2 and a. 


The Mach number in Region 2 can be solved for readily since the flow is assumed to be 
isentropic during expansion. Total (stagnation) pressure in the nozzle is measured with a 
Pitot tube on the thrust stand located at the entrance to the nozzle. In the isentropic 
expansion of an ideal gas, Mach number and total to static pressure ratio are related as 
follows. 


P/p = (1.0 + 



a-nn 

M 2 ) 


Eqn(10) 


In Equation 10, (P) is the total pressure, (p) is the static pressure, (M) is the Mach 
number. Solving for the Mach number. 


2 .A 

M = ( [ (P/p)<-H)/7 . 1.0] )' A Eqn(lOa) 

7-1 

Since the final Mach number is known, Equation 1 can be solved for the initial flow 
direction of the jet boundary. Initial boundary flow angle increases as inlet pressure 


increases. 


47 


Principles of Expansion Wave/Oblique Shock Wave Interaction 

The most important concepts associated with expansion wave and oblique shock 
wave reflection and interaction are as follows Ref(4). 

Expansion Waves: 

1. An expansion wave turns the flow away from itself. 

2. Expansion waves reflect from solid surfaces as expansion waves. 

3. Expansion waves reflect from a free pressure boundary as oblique shocks. 

Oblique Shock Waves: 

1. An oblique shock wave turns the flow toward itself. 

2. Oblique shock waves reflect from a solid surface as oblique shocks. 

3. Oblique shocks reflect from a free pressure boundary as expansion waves. 

Derivation of the Maximum Flow Deflection Angle (g m av) 

The Mach Disc forms when an oblique shock wave with a large deflection angle (5) 
generates a reflected shock with a corresponding deflection angle greater than (S)max- 
When the deflection angle corresponding to the Mach number in Region 5 (behind the 
incident oblique shock) exceeds (5),^, the Mach Disc will form. The maximum deflection 
angle associated with an oblique shock is obtained in the following manner. The geometry 
of flow through an oblique shock wave is illustrated in Figure 9. 


48 



Figure 9. Flow Geometry Across an Oblique Shock Wave 

Upon inspection of Figure 9, the flow deflection angle (£) can be defined in terms of 
velocities as follows. 

5 = tan- 1 (V,/V v) - tan-' (V t /Vn,) 


Where V t and V n are tangent and normal velocity components, respectively. Recall for an 
oblique shock wave that the tangential velocity remains constant across the shock. (g maT ) 
can be calculated by differentiating ( s ) with respect to V t in the above expression and 
setting it equal to zero. After differentiating, 6,^ occurs when: 

V t /V n| = CVn 2 /V„,)^ 

Substituting into the equation for 6. 

<5 max = tan -1 (Vni/Vn 2 )^ - tan- 1 (Vni/Vn 2 ) 

From the governing equations of oblique shock waves. 

Vni/Vn 2 = pi! p\ = r 


Substituting 




^ max = tan ' 1 (r)' / * - tan- 1 (r)' l/4 
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Employing a trig identity and rearranging. 

5 max = tan -1 - — — Eq n (5) 

2(r) 2 

If the flow turning angle associated with a given Mach number is greater than 5,^, then a 
regular Mach reflection is not possible and a Mach Disc consequently forms. Ref(7) 
presents an excellent development of Mach Disc or Mach reflection. 
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Appendix B 


sample Calculation (No Mach Disc) 
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Sample Calculation of Jet Structure (No Mac h Disc) 

In reference to Figure 2, the shock structure within the constant pressure boundary 
of an underexpanded jet can now be calculated. In this sample calculation the ratio of 
nozzle exit static pressure to ambient pressure is 2.0. For a constant specific heat ratio of 
7 = 1.4, Equation 3a can be written as follows: 

M 2 = [ (5.0((P/p) 0 - 2857 - 1.0) ]' A 


For pj/p 0 = 2.0, pj is then 29.4 psia. In order for the flow to choke, the pressure must 
exceed the critical pressure ratio of P/p = 0.52828. The flow will therefore be choked for 
any nozzle total pressure greater than 27.8 psia which is the case in this example. For 
p./p o = 2.0, P/p - 3.8. Solving Equation 10a for M 2 . 


Substituting this value of Mach number into Equation 2, the initial turning angle of the 
flow is 12.5°. The next step is to determine the angle of the initial Mach wave. 



Figure 10. Prandtl-Meyer Expansion 
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Recall from Prandtl-Meyer expansion Ref(2). 


tan( M ) = 1.0/ (M 2 - 1.0) 1/4 
or 

m = tan 1 [ 1.0/CM 2 - 1.0)* ] Eqn(ll) 

In Equation 11, 0) is the Mach angle. For Mi = 1.0 using Equation 11, (mi) = 90° and 
M 2 = 1.52, (m 2 ) =41.1°. From Figure 10, the average of (mi) and (m 2 ) measured from 
the horizontal axis is: 


Mave 


mi + (n - Q ) 
2 


Eqn(12) 


Substituting values (M ave ) = 59.3°. This is the angle the first expansion wave makes with 
respect to the horizontal axis. Referring once again to Figure 2 as the flow travels from 
Region 3 to 4 it passes through a reflected expansion wave. From Prandtl-Meyer 
expansion theory the flow deflection angle (a) is related to the initial and final Prandtl- 
Meyer angles as follows. 


a = 1/2 - ^1 


or this case: 


«24 = V4 - V2 


(i/ 2 ) was calculated before from Equation 1 as 12.5°. The initial flow turning angle was 
12.5° so (1/4) was calculated as: 

i/ 4 = 12.5 + 12.5 = 25°. 
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with (i/ 4 ), {v 2 ), and (M 2 ) known, Equation 1 can be solved implicitly for M 4 . 

M 4 = 1.95. 

The average Mach angle or discrete Mach wave for Region 2 to 4 was solved in a 
similar as that from 1 to 2. However, it is important to keep track of the relationship 
between these angles and the flow axis. 



Figure 11. Expansion From Flow Region 2 to 4 


Substituting values for M 2 an M* into Equation 11, («) = 41.1 + 12.5 = 53.6°, Om) - 
43.2°. The average Mach wave from Region 2 to 4 is 48.6° with respect to the horizontal 
axis. Recalling the shock structure from Region 2 to 4. 



V 


Figure 12. Jet Structure From Flow Region 2 to 4 
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In Figure 12, the non-dimensional distance xi/r„ is easily solved for as xi — r„/tan(59.3) 
or xi/r n = 0.593. The distance \i!x n can be solved for, knowing the equation for the 
curved boundary. The jet boundary solution for pj/p 0 = 2.0 from the Jet Boundary 
Fortran program was put into a curve Fitting program and the following 3rd order 
polynomial was generated. 

r/r n = 1.0 + 0.225(x/r n ) - 0. + 0.033(x/r n )’ Eqn(13) 

There are 2 unknowns in Equation 13, r/r n and x/r n . Another expression relating x/r n and 
r/r n can be developed upon inspection of Figure 12. 

x/r n = xi/r n + x 2 /r n = 0.593 + x 2 /r n 
x 2 /r n = (r/r n )/tan(48.5) 

Combining both equations to get one equation in terms of r/r n and x/r n : 


r/r n = tan(48.5)(x/r n - 0.593) Eqn(14) 

Combining Equations 13 and 14 and solving for x/r n and r/r n : 
x/r n = 1.61; r/r n = 1.15 

Up to this point the flow is assumed to be isentropic since it has only encountered 
expansion waves. The static pressure at location 4 can then be calculated using isentropic 
and expansion process. 

(Tl)/(1) 


P4/P4 = (1.0 + ( (7 - l)/2 ) M4 ) 


Eqn(15) 
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Where (7) = 1.4, M* = 1.95, P 4 = P n = 55.65 psia. Substituting values (p 4 = 
7.69 psia). Recall from the general principles of expansion waves that an expansion wave 
reflects from a free pressure boundary as an oblique shock wave. The static pressure in 
Region 4 is less than atmospheric. The flow must therefore encounter an oblique shock 
wave in order to raise its pressure to atmospheric pressure in Region 6 and maintain the 
constant pressure boundary. 



Figure 13. Oblique Shock Wave Structure 


In Figure 13, (5) is the flow turning angle, (e) is the oblique shock wave angle with 
respect to the original flow axis, and (fi) is the angle of the oblique shock wave with 
respect to the flow direction on the other side of the shock wave. From Ref(2), the 
governing equations for property change across oblique shocks are as follows. 

27 7 - 1 

p 6 /p 4 = M 4 2 sin 2 ( £46 ) Eqn(16) 

7+1 7 + 1 

Static density ratio: 


(7 + 1) M 4 2 sin 2 ( £46 ) 


P6b* — 


2 + (7 - 1) Ma 2 sin 2 ( £4 6> 


Eqn(17) 
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Relationship between flow turning angle ( 5 ) and oblique shock angle ( c ): 

(7 + 1 ) W 

l/tan(6) = [ - 1 ] tan(c46) Eqn(18) 

2 M4 2 sin 2 ( C46 ) - 1 

Mach number on the other side of an oblique shock: 

tan(c46) 2 1 7 - 1 

= ( ) ( r + ) Eqn(19) 

tan 09 4e) 7+1 M$ 2 sin 2 (pu , ) 2 

In the above equation, the unknowns are (c46 ), (p6), ( 646), (Me) and (£46 )• In most 
problems involving oblique shock waves there is something known on the other side of the 
shock. It is usually the flow turning angle (6). In this case the static pressure in Region 6 
is known and is atmospheric. The procedure for solving for the unknowns from the above 
equations is as follows. 

1. Solve Equation 16 for the oblique shock wave angle (c46)- 

2. Solve Equation 18 implicitly for the flow turning angle (546)- Solve for 
(fix) = (£46) - (546). 

3. Solve Equation 19 implicitly for M6. 

For (7 = 1.4) Equation 16 can be rewritten as follows. 

P6/p4 = 1.1666 M4 2 sin 2 (£46) - 0.1666 


Where p$ = 14.7 psia, p4 = 7.69 psia and M4 = 1.95. Solving for (c «), (c46) = 43.2° 
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Solving Equation 18 for (546) and using (7 = 1.4). 

1.2 M 4 2 

546 = tan 1 [ ( - 1 ) tan(c4e) ] * 1 

M 4 2 sin 2 (c46) - 1 

(546) = 12-42°; 0?46) = (£46) - (546) = 30.8° 

Solving Equation 19 for M«. 

tan (£46) 

M 6 = [ ( - 0.2) sin 2 ( 046 ) ] = 1.50 

0.833 sin 2 (fi 46 ) 

In this region of the jet, two oblique shock waves of equal strength intersect at the 
centerline of the flow stream. When oblique shocks of equal strength intersect, the 
problem is treated in the same manner as a normal reflection from a solid surface Ref(l). 
The flow deflection angle (567) cannot exceed the maximum flow deflection angle (5) max . 
Recall from the discussion of Mach Disc formation and Equation 5 that (5)5^ is defined as 
follows. 


r - 1 

5™, = tan- 1 r = P 1 I P6 

2(r) 2 


and (p7/p6) from Equation 17 is 


(7 + 1)M 6 2 sin 2 (£67) 

T67 = Pi! P6 = 

2 + (7 - 1) M 6 2 sin 2 ( £ 67) 

In this equation there are two unknowns (£67) and (^7). (£67) corresponds to the maximum 
oblique shock wave angle (e),^ and is found by differentiating Equation 18, setting it 
equal zero and solving for (e),^. 
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sin 2 ( £max ) 


1 (( 7 + 1) M^-l) [ (7 + 1) (1 + (7 - 1) M6 2 7 + 1 

{ + + } 

7 M* 2 4 2 16 

Eqn(20) 


Substituting for (e),^ = 66.5°. Solving Equation 12 for r 67 and solving 
Equation 5 for ($)max- 

r 67 = 1.617 (i) m = 14.5° 

The value of (6)^ is greater then (5 67 ), 12.5° and consequently there is no 
formation of a Mach Disc. The reflected shock from Region 6 to 7 is therefore strong 
enough to straighten the flow parallel to the centerline axis. The only known property in 
Region 7 is therefore the flow deflection angle. That is ( 545 ) = (56 7 ) = 12.42°. Looking 
at Equations 16 through 19 and knowing (6 67 ), Equation 18 can be solved implicitly for 
(cg 7 ). ( C67 ) = 65.3°. Substituting into Equation 16 and solving for pi. 


pi — P6 (1.166 M^ 2 sin 2 (cei) - 0.166) 
P "7 = 29.5 psia 


Note that p 1 is higher than atmospheric pressure which is expected since the flow passed 
through an oblique shock from Region 6 to 7. The final step in this problem is to calculate 
where the reflected shock wave will intersect the jet boundary. Figure 2 is repeated at this 


time for clarification. 
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Figure 2. Jet Structure of an Axisymmetric Jet Exhausting Into Still Air 

Recall from the previous discussion that the inflection point occurs where the 
expansion wave reflects from the free pressure boundary. Beyond this point, the curvature 
of jet boundary is assumed to be the mirror image of the boundary up to the point of 
inflection. Recall Equation 13. 

r/r n = 1.0 + 0.225 (x/r n ) - 0.137 (x/r n ) 2 + 0.033 (x/r^ 

The non-dimensional distance (xi + x 2 ) was calculated earlier as (1.61). Rewriting 
Equation 13 so that the boundary beyond the point of inflection is the mirror image of the 
boundary before the inflection point. 

r/r n = 1.0 + 0.225(1.61 - x/r n ) - 0.137(1.61 - x/rj 2 + 0.033 (1.61 - x/r n ) 3 Eqn(21) 


where 


x/r n = x 3 + X 4 
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The value of r/r n is the value at the inflection point and was calculated earlier 

(r/r n = 1.355). 

From Figure 12: 

X 3 = (r/r a )/tan(c 46 ) Eqn(22) 

X 4 = (r/r^/tan^a?) Eqn(23) 

Where 067 = £67 - 567 = 65.3 - 12.5 = 52.9 

Substituting X 3 into the x/r n expression 

x/r n = X 3 + X 4 = (r/r n )/tan(c 46 ) + X 4 Eqn(24) 

The value of r/r n at location 3 is the value at the point of inflection (r/r n = 1. 145). Recall 
that (c 46 = 42.3°) 

Solving Equation 22 for X 3 

x 3 = (1.145)/tan(43.2) = 1.219 

Substituting X 3 into Equation 24 

x/r n = 1.219 + x 4 Eqn(25) 

Looking at Equations 21, 23 and 25 there are 3 unknowns, x/r n , X 4 and r/r n . 

Solving these equations simultaneously for the unknowns. 

X 4 = (0.335) r/r n = (1.069) 

x/r n = x 3 + x 4 = 0.335 + 1.219 = 1.554 


x/r n (total non-dimensional distance in x direction) = xi + x 3 + X 4 


x/r n (total) = w/r n = 1.61 + 1.554 - 3.16 


The distance (w/r n ) represents the point at which the reflected oblique shock wave 


meets the free pressure boundary. 


Appendix C 


Sample Calculation (Mach Disc Present) 
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Figure 4. Mach Disc Structure 


In Figure 4b, Regions 4, 5, and 6 coincide those in Figure 2. The values of static 
pressure, temperature and Mach number in Regions 4, 5, and 6 are calculated in the same 
manner as before. As the flow travels from Region 4 to 9 it passes through a normal 
shock. The following are the governing equations for property change across a normal 
shock wave. 


Static pressure ratio: 

( 2 7 ) W 

P9/P4 = -j— j 

7 + 1 



Eqn(26) 


Static density ratio: 


P*)! P4 = 


(7 + 1 ) MS 
2 + (7 - 1)M4 2 


Eqn(27) 
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Static temperature ratio: 


Wu — 


2 7 M4 2 



7-1 
7 + 1 


)( 


7-1 

7 + 1 


+ 


2 

( 7 + 1 ) M4 2 


) 


Eqn(28) 


Mach number on the other side of a Normal Shock: 


M4 2 + 2 /( 7 - 1 ) 

M 9 = [ ] 

(2 7 / ( 7 - 1 )) M4 2 - 1 


Eqn(29) 


As seen in Figure 4b, the Mach numbers in Regions 7, 8, and 9 are different. The 
only flow property in common between Regions 7, 8, and 9 is the static pressure. From 
symmetry the flow properties in Regions 7 and 8 will be the same. The procedure for 
calculating the flow properties in Regions 7, 8, and 9 is as follows. 

1. Solve for flow properties in Regions 4, 5 and 6 as before (ps = P4 = p a tm)- 

2. With p4 and M4 known, use Equation 26 to solve for ps. 

3. Use Equations 27, 28 and 29 to solve for p 9 , t 9 and M 9 . 

4. Since p7 = p 9 use Equation 16 to solve for the oblique shock wave angle (£57). 

5. With M 5 and (c 57 ) known, use Equation 18 to solve for the flow turning angle (557) 
implicitly. Verify that this angle is less that (S)^^ at that nozzle pressure. 

6. Use Equation 19 to solve for M7 implicitly. Recall that ifisi) — (^57) - (657). 

Static pressure, temperature, density and Mach number in all regions up to Region 
9 are now known. The next step is to determine the size and location of the Mach Disc. 
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w 


If the Prandtl primary wavelength (w/r n ) is known by actual measurement from a 
Schlieren photographs or by use of Equation 4, then the structure of the Mach Disc can be 
determined using the data calculated in the Jet Boundary and Jet Structure Fortran 
programs. The flow properties in Regions 4, 5, 7 and 9 were calculated usmg the Jet 
Structure program. From symmetry, properties in Regions 6 and 8 are equal to those in 5 
and 7, respectively. The non-dimensional distance x/r n is also known for the case of 
nozzle total pressure of 89.0 psia (pj/ Po = 3.2). The calculated properties are as follows. 


(£57) = 42.4 c 


(£45) = 43. 2 C 


x/r„ = 1.99 


w/r n = 4.05 (Eqn 4) 

w/r n = 3.95 (measured from Schlieren Photo No. 30) 

If the Prandtl primary wavelength is taken as 4.05 then from Figure 4b. 
X4 + X5 = w/r n - x/r n = 4.05 - 1.99 = 2.0 


X 4 = r4/tan(£45) 


X 5 = rs/tan(£57) 


rm4 = r 4 + r d = 1.355 (from Jet Structure Program) 


rms = rs + r d 


Eqn(30) 

Eqn(31) 

Eqn(32) 

Eqn(33) 

Eqn(34) 


There are 6 unknowns and 5 equations. The unknowns are X4, X5, r4, rms and 
r d . However, since the angles £45 and £57 are approximately equal the non-dimensional 
radii rrm and Tms will be approximately equal. 

T m4 = F m5 = 1.355 
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If rn»5 is known then Equations 30 through 34 reduce to 5 equations and 5 unknowns. 

The procedure for solving these equations is as follows. 

1. Substitute Equation 31 and 32 into 30. 

2. Solve Equations 33 and 34 for r 4 and and substitute into Equation 30. 

3. Solve new Equation 30 for r d . 

4. Solve for the remaining unknowns. 

The values of the unknowns are as follows. 

X4 = 1.017 x 5 = 1.047 

r 4 = 0.955 r 5 = 0.955 

and r d = 0.4 

The value x ra + xa is the distance from the nozzle to the Mach Disc. In this case: 

L d = x m + *4 = + = 3-007 

This value is also measured from the Schlieren photograph and is 3.0. This also agrees 
with the NASA data in Ref(13). 


Appendix D 


Schlieren Photographs 


67 


Schlieren Photographs 


P n = Nozzle Total Pressure 
pj = Nozzle Exit Static Pressure 
p 0 = Atmospheric Pressure (14.7 psia) 


Photo No. 

P n (psia) 

P/Po 

1 

32.4 

1.164 

2 

35.4 

1.272 

3 

41.9 

1.506 

4 

51.8 

1.862 

5 

57.8 

2.077 

6 

66.2 

2.379 

7 

67.5 

2.426 

8 

75.2 

2.703 

9 

79.8 

2.868 

10 

88.7 

3.188 

11 

96.1 

3.454 

12 

32.0 

1.150 

13 

35.0 

1.258 

14 

42.1 

1.513 

15 

53.0 

1.905 

16 

58.7 

2.109 

17 

66.3 

2.383 

18 

68.0 

2.444 

19 

76.6 

2.753 

20 

79.8 

2.868 

21 

89.5 

3.216 

22 

31.9 

1.146 

23 

35.6 

1.279 

24 

41.2 

1.481 

25 

52.2 

1.876 

26 

57.8 

2.077 

27 

66.5 

2.389 

28 

68.4 

2.458 

29 

75.6 

2.717 

30 

88.6 

3.184 

31 

97.0 

3.486 

32 

32.1 

1.154 

33 

35.5 

1.276 

34 

41.7 

1.499 

35 

52.0 

1.869 

36 

55.2 

1.984 





Schlieren Photo No. 18 
= 68.0psia 


Pn 







Schlieren Photo No. 19 


?n = 76.6 


psia 
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Calculated Jet Boundary Contours 


















Appendix F 


Jet Boundary Fortran Program 


on 




PROGRAM 1 




r >:/- 


\_ “"T- s. 


PROGRAM 

INTEGER 


PR0G1 

I 


REAL A f B ( 420 ) , DR < 420 ) , AST , ALF,TH,R,RR, AE , BRO , BRA , BRB , 

BRD , FRO , FRA , FRB , FRC.FRD , MRA , MRB.MRC, HRB , AA , MACH , 
CR<420>,BA<420>,DD,DF,PTJK,PTJ,BME,PJH,PJ,ME,BRC,H 
DATA B/420*0.0/ - 

DATA DR/420*0.0/ 

DATA CR/4 20*0.0/ 

DATA BA/420*0 . 0/ 

PRINT*,/ ENTER NOZZLE 'TOTAL PRESSURE* IN PSIA' 
PRINT*, * ENTER 0.0 IF COMPLETED* 

READ*, PTJ 
IF (PTJ.EQ.O.O) THEN 
GOTO 300 
END IF 

H = 0.00108333 
PJ = PTJ*0. 52828 
PJF: = FJ/14.7 
F'TJR = PTJ/14.7 

ME = SQRT(5.0*PTJR**0. 285714-5.0) 

BME = SORT (ME**2-1 .0) 

ALF = 2. 449489*ATAN<BME/2. 449489) - ATAN(BME) 

AE = < 1 . 0/ME > * < 0 . 8333333+0 . 1 6A66£6*KF.**2 ) *#3 
R = 0,2275 + H 
FRO = 1.0/TAN (ALF) 

DR(1) « 1.0 

DD =0.0 

DF = 0.0 

CR ( 1 ) = 0.2275 

PRO = SORT ( ME**2-1 .0 ) 

PRINT*, 'ASTE = * , AE 


DO 10 


PRINT*, 'ME = '.ME 
1 = 1 .35 

RR = R/0.2275 
AST = (RR**2)*AL 
MRA = MACH (AST) 

PR A = SORT ( MRA**2-1 ♦ 0 > 

TH = 2 . 449489* ( AT AN < BRA/2 . 44948? > -ATAN ( PRO/2 . 

ATAN(PRA) + ATAN(BRO) 

FRA = 1 .0/TAN<ALF-TH) 

F: = R + H 
RR = R/0.2275 
AST = <RR**2)*AE 
MRB = MACH (AST) 

BRD = SORT (hRB**2-l .0) 

TH = 2.449489*(ATAN( BRB/2 ♦ 449489) - A FAN ( &F;0/2 . 

ATAN (BRD) + ATAN (PRO) 

FRB = 1.0/TAN(ALF-TH) 

R = R + H 
RR = R/0.2275 
AST = <RR**2)*AE 
MRC = MACH (AST) 

BRC = SORT <MRC**2— 1 .0) 

TH = 2 . 449439* (ATAN < BRC/2 . 449489 ) -ATAN ■ PRO/2 , 
ATAN (BRC) + ATAN (PRO) 
i-RC = 1. 0/TAN (ALF-TH) 

A = 0 . 375*H* ( FRO 1 3 . 0*FR A4 3 . 0*FRB+FRC ) 

AA = A/0.2275 
DR(T-M) = RR 
DD = DD + AA 
B ( J -! 1 ) = PD 
CR( 3+1 ) = R 
DF = DF + A 
BA ( 1+1 ) = DF 
FRO - r RIi 
R = R 1 H 
CONTI Nur. 


449489) )- 


449489) )■ 


449409) )- 




o* 


0 

9 


O* oO' 


A® 






10 


Wf-J 


> 'WB& Boiiti WrfBIS’ft ASiiEtf ’ Assure. W «o^pjr 

•:> ~ : - ‘X -?. '■** ~~v- S- vV; •' ' X^trZr- ■ ' • 

PRINT*, J W fe- V ,H'"-: — • -- ;• • -; '...'. :.. .. 


SOLUTION FOR JET BOUNDARY' 
R/RN X 


PRINT* 

PRINT*,' 

PRINT*^' R 

D ° *PRINT21 » 'CR< I ) ,DR< I ) »BA< I ) ,.B< I ) _ 

1 FORMAT <ix,F6.4,3X,Fl2.9.3X,F12.8,3X,F12.8) 

0 CONTINUE 

GOTO 5 

300 REAL FUNCTION MACH<xV w 

REAL MU , ML , HA , FMU , FML , FMA , X 

100 FML <l«0/ML)*<0»8333333+0. 1666666*ML**2)**3 - X 

FMU = = M ( 1 !o/MU?*(O.B333333+0.1666666*MU**2)**3 - X 

FMV = FMU*FML 
IF<FMV.GT.O.O> THEN — 

ML = MU 
GOTO 100 

END IF 

FMA = **<1 .0/MaT*U 0^8333333+0. 1666666*MA**2)#*3 - X 

FMV = FML*FMA 

IF< ABS(FMA) .LT.0.1E-4) THEN 
MACH = MA 

£pyj |n^| 

ELSE IF <+ mO.LT. 0.0) THEN 
MU = MA 

F MU = ( 1 .O^MU)* < 0 . 833333+0 . i 66££»66*MU**2 ) **3 - X 
COTO 200 

ELSEIF (FMV.GT .0.0) THEN 

FML~= M ('l .0/ML)*(0. 833333310. 1666666*ML**2)**3 - X 

GOTO 200 

END IF 
END 


X/RN 


END OF FILE REACHED - PRESS RETURN 
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Jet Structure Fortran Program 


PROGRAM 2' ' ,r ^ ’ 
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C 

c 

c 

c 

c 

c 

c 

c 


PROGRAM 

REAL 


PR0G5 

PT J f PT JR r P J r MA , MUW T BME , BEL , ALF , XA , NU , ML. , MU , FML , FMU , 
FMV MB , FMB , HUB , MUAV , FXL , XL , XU , EX , FXA,FXV , FXU , EBC , PB , 
RRN , A , B , C , D , XC , DEL B , BET , MD , XD , MUA , BMU , BMB , EL , t.A , EU , 
FEA,F£u,FEV,FEL r ECP,FXCV,FXCU,FXCL,FXC,XRN,BML,EM,XE, 

BETB,RPR, AB 

TOTAL PRESSURE IN PSIA' 

COMPLETED' 


dl 1 j r\jur> » nr 

PRINT*, 'ENTER NOZZLE 
PRINT*, 'ENTER 0.0 IF 
READ*, PT J 
IF <PTJ . EQ .0.0) THEN 
GOTO 950 
ENBTF 

PRINT*, ENTER 
READ*, TN 
PRINf *, 'ENTER 
READ*, A 
PRINT*, 'ENTER 
RFAD*,B 
PRINT*, 'ENTER 
REAP*, C 
PRINT*, 'ENTER 
READ* , P 

PJ = PTJ*0. 52828 
F* JR = PJ/14.7 
PT JR = PTJ/14.7 

MA = SORT ( 5. 0*PTJR**0. 285714-5.0) 

MU = ATAN(1.0/( <MA**2. 0-1.0) **0.5) 
BME -• SOF T ( MA**2 . 0- 1.0) 

_ \T AN ( 


NOZZLE TOTAL TEMPERATURE IN DEG. RANKING 

COEFFICIENT A' 

B' 

C' 

B' 


COEFFICIENT 

COEFFICIENT 

COEFFICIENT 


DEL. = 5 . 44949?*ATAN ( BME/2 . 449489 ) -A I AN ( BHC ) 


ALF - < 


570796+MU DED/2.0 


XA - 1.0/TAN (ALF) 
MUW = AT AN ( (MA**2 • 


111 A 


MUW + DEL 


> - ■ 1 . 0 ) S* ( -o .5)) 


PRT 

Ir t y 

7 Mfi = 

' ,MA 

PR I 

I'll ^ y 

'PEL = 

' , DEL 

PR I 

NT*. 

'ALF = 

' , ALF 

pr: ; 

NT# * 

'XA = 

' f XA 

pf;i 

NT# r 

'h-JW - 

= 7 , MUW 

PR I 

NT# f 

'MUA = 

= ' f MUA 


NU = 2 . 0*DEL 
ML =1.0 
lr ML = SORT (ML**2 .0-1.0) 

FML = 2. 449489*ATAN (BML/2. 449489)- 


-AT AN ( BML ) -NIJ 


MU = ML 1 0.25 


BMLS = SOFT (MU**2. 0-1.0) 

FMU = 2. 449409*ATAN(BMU/2. 449489) 

FMV = F MU* FML 

IF < FMV . GT .0.0) I HEN 

ML = MU 

GOTO 100 

END IF 


-A J AN • BMiJ ) -NU 


MB = (ML 

BMB = SORT (MB**?. 0-1.0) 

FMB = 2.4 4?48?*ATAN ( BMB 72 . 449489 ) -AT AN ( BMB ) - NU 
FMV = FML*FMD 

IF ( ABS ( FMB) .LT.0.1E-4) THEN 
PRINT*, 'MB = ' .MB 
GOTO 300 

EL.SEIF ( MV . LT . 0 . 0 


HEN 


MU = MB 


GOTO 200 

ELSEIF ( PMU . GT , 0 . 0 ) THEN 
Ml = MB 

BML - SORT (,"'..**2 . 0-1 . •“•> 
FML = 2 . *• 4 r 4 8 9 * A '( A N ■ 
GOTO 200 
END' F 






/!49<R9 


! AN ( ?•" 










400 


500 


600 


800 


PRINT * . ' NUB - , NUB 

PRINT*; 'MUAV === SMUAV 
XL = 0.0 

FXL = A-TAN ( MU AO ) * < XL-XA ) +B*XL ! C*XL**2.0+D*XL**3 . 0 

pXU~= X k-TAN ( MUAV ) *(XU-XA> +B*XU ! C*XU**2 . 0+D*XU**3 . 0 

FXV •'* FXU*FXL . ;. ' . . . . . - 

iF(FXV.GT.O.O) THEN 
XL = XU 
GOTO 400 
END IF 

XD = (XU+XD/2.0 

FXD = A-TAN ( MUAV ) * ( XD-XA > +B*XD+C*XD**2 ♦ 0+D*XH**3 . 0 

I^yABS^FXD^llT.O.lE-4) THEN 
PRINT*. 'XD = ' »XD 
GOTO 600 

I^SEI^FXV.LT.O.O) THEN 
GOTO 500 

ELSEIF(FXV.GT.O.O) THEN 
Xf = XD 

FXL = A— TAN(MUAV)*(XL-XA) +B*XL+C#XL**2 .04 D*XL**3.0 

GOTO 500 

ENDIF 

IIcIaII $ 1 5 tK?*38B8?i2J84i?*? 1 . o/ < i . » ^^m***;- . o > > . «« . 

PRINT*, 'EEC = ' , EBC 

PRINT*. 'PB = ' , PB 

RRN = A+B*XD ! C*XD**2 . 0+D*XD**3 . 0 

XC = RRN/TAN<EBC> 

PRINT*, 'RRN = ' ,RRN 
pp. j jj/ x XC = y f XC 

DFl P = J AT AN ( ( (1 . 2*MB**2 • 0/ ( ?1B**2 . O* (SIN ( EEC ) ) **2.0-1 .0) 1 .0)* 
TAN(EBC) >**<-1 .0) ) 

PET = EEC - DELE 
PRINT* , ' DELB = ' , DELB 

Vili 1 = 7 ( l TAN ( FBC ) >(b ; 8333333 *T AN ( BET ) >- 0 . 2 ) * ( SIN < BFi ) )**2.0) 
**(-0.5) 

PRINT*, 'MD = ' , MD 

PRINT* 

DM - (<1.7117)*PTJ)/TN 
DN*1. 577440? 

DTN/ ((1.0+0. 2*MB**2 . 0 ) **2 . 5 ) 

( ( MB**2 . 0+5 . 0 ) / ( 7 . C*MI<**2 . 0- 1 . 0 > ) **0 . 5 
PB* (1,1 6 6 666 6 * H B * * 2 . 0-0 . 1 66666/ ) 

PF 

DB* ( 2 . 4*MP**2 . 0* ( ST N < t PC • ) **2 .0/(2. 0+0 . 4*KB**2 . 0* 

(SIN (EEC) >**2.0) 

3? .6? 17/ED 

DB* ( 2 . 4*MB**2 . 0 ) / ( 2 . 0+0 . 4*MB**2 . 0 ) 

< 2 ♦ 60251 ?*PF ) /DP 
= ( 0 . 086525*PTJ) / ( TN**0 . 5 ) 

' ,TN 
' r I«N 
' , DB 
' , r«F 
' t PF 
' , DD 
STD 
' , DF 

= ' , HDOT 

EM - OSIN < ( < 1 .0/(1 . 4*MD**2. 0) ' * ( 0 . 6*MD**2 . 0-1 .01(2.44 
0 . 48*MD**2 . 0+0 . 36*MD**4 . 0 ) **0 . 5 ) ) **0 . 5 ) 

ERA = <2. 4 *MD **2 . 0# ( SIN (FH> > **2 . 0 > / <2 . 04 0 . 4*MD**2 . 0* 

( SIN (EM) >**2.0) 

UMAX - ( DRA-1 . 0 ) /( 2 . 0* ( DRA**0 . 5 ) ) 

PR 1 N f * . ' EMAX= ' , EM 
PRINT*, 'DRA = ' , DRA 
PRINT*, 'UMAX = SDMAX 
IF ( DELE. GT. DM AX) THEN 
PR INI*, ' **** *DMAX IS EXCEEDED*****' 

UR J - 2 . 0* ( 1 . 52* ( ( F JR ) *#0 . 4 3 7 ) -0 . 5 ) „ ... 

PRINT*, 'PRANDTl WAVELENGTH (NASA) j 5 
ELSE 

PRINT*. '*****DMAX IS NOT EXCEEDED*****' 

UR.: =• 2 . 0* ( i . 55* ( ( P Jfi 1.0) **0 . 5 > ) 

... ' , I •lAIIVI l-M.-.TII 1 I * / A \ 


DTN = 

DB = 

HF = 

pp — 

PF = 

DD = 

TD = 

DF - 
TF = 

HDOT 

PRINT*, 'TN 
PRINT*, ' DN 
PRINT*, ' DB 
PRINT*, ' HF 
PRINT* * 'PF 
PRINT* , ' DD 
PRINT*, 'TD 
PRINT*, 'DF = 
PRINT* , 'HDOT 


PRINT * f 
FNDIF 


PRANDTL yAVUL CNGTH < Nr.Or^ ) 


. UR J 




o< 


w 


. . . . . . r- f , I f > ' 

. i iL - f*.--!;-* I • v }?(;♦ ] 66^001 •// v 1 .« : L* • \, ; ; : « V-V* 4 L> ' 

DE = ( DD#2 . 4*MD**2 . 0* (SIN ( EDE ) ) **2 . 0 ) / ( 2 . 0+0 . 4*MD**2 . 0* 

C HELD = S ATAN? < < < 1?2*«D**2 . 0 ) / ( MD**2 . 0* ( SIN ( EDE ) ) **2 .0-1,0) 

C -1 .0)*TAN(EDE) )**(~1 .0) ) 

BET * EDE - DELD 
ME *' < ‘ 


**<-0.5) 


DE ) / CO . 8333333*T AN < BET) > -0 . 2 ) * < S I N < BET ) ) **2 .0 ) 


I F < DELB . GT . DMAX ) THEN 
PRINT* 

RH = 1 . 0/T AN < EBC ) + 1 . 0/TAN < BET) 
RF = < 1 ,0/RH)*(WRJ-XD> 

PRINT*, 'RF = ' ,RF 
PRINT*, 'RH = ' ,RH 


TE = <2.70012*PE)/DE 
TB = <2.70012*PB)/DB 
TF = <2.70012*PF)/DF 


C 


RD = 
XX = 
XZ = 
LD = 
XRN = 
LDN = 
RDN = 


RRN - RF 
RF/TAN(EBC) 

RF/TAN(BET) 

XD + XX 
vr, i. yy 4- Y7 

0 . 055038+1 . 396*P JR-0 . 1978* < P JR**2 . 0 ) -I 0 . 0 1 3 1 4* < P Jl .**3 . 0 ) 
- 0 . <?4807 +0 . 67296*P JR-0 , 08638* < P JR**2 . 0 ) +0 . 004326* 
<PJR**3.0) 


l HS=( DB*MB*49 . 032* ( TB#*0 . 5 ) *3 . 14159265* ( ( RRN*0 . 2275 ) **2 . 0 ) ) / 1 44 < 
RZ = (TAN (EBC) )*XX 
RDB = RRN - RZ 

AB = 3.141 592654* < ( RRN*0 . 2275 ) **? . 0 > 

2 .3 41 592654* < < RI«B*0 . 2275 ) **2 . 0 > 

ABC - AB/ < S3 N < EBC ) ) 

fV - 3 . i4 1592654* < < RPB*0 . 2275 ) **2 . 0 ) 

* E TB _ pjxp — rjjjj p 

RHS = < DD*MD*49 .032* ( TD**0 . 5 ) * (SIN ( BETB ) > *ABC ) / 144 . 0 
1 < RF *MF *49 . 032* < TF**C . 5 ) *AF > / 1 44 . 0 


PRINT*, 

I R I ! '■’ • 'J n 

PRINT*, 
PR I NT * , 
PRINT*, 
PRINT*, 


f 7. = '.XZ 
'DISTANCE TO 


DISC IF? 


DISC (NASA 

is: ' r RD 
1 <? ♦ 


Jl'rf. J ! IIHWL 1 W MACH 
DISTANCE TO MAC!! 

MACH DISC RAD i US x ^ 

KKJ.Nl*, 'MACH DISC RADIUS ( NASA ) IS t 
PRINT*, 'CALCULATED PRANDTL PRIMARY 
PERC = (LDN/LD )*100.0 
PRINT*, 'PERC (LD) = ' 

PRINT * , ' UR J = ' ,WRJ 

MACH 


: ‘ERC 



PRINT*, 'PERC (RIO = ',PERC 
ERR = APS ( 100 . 0— PERC ) 

PRINT*, 'PERCENT DEVIATION OF MACH 
PRINT* 

, LHS 


I JtltT I <?• 

PRINT*, 'LHS = 

PRINT*. 'RHS = 

PRINT*, 'RZ = 

PRINT*, 'AB = 

PRINT*, 'RDB = 

PRINT*, 'ABC = 

PRINT*, 'AF - 
END IF 

PRINT*, ' BET2 
PRINT*. 'ME = 

PRINT*, 'TE = 

PRINT*. 'T # = 

PRINT*; 'TF = , T 

IF (PELB.Cn DM AX) 

GOTO 60 

FNBT F 

RRNB - A + b * ( XD-XC ) 3 C* ( XD-XC > **2 - V+I<> < XD-XC ) *4 3 . 0 
XF - RPNB/I AN 'PDF) 

-- >:c + xe 


, i-D 

: is: DM 

' , RDN 

WAVELENGTH IS: '.XRN 

DISC DISTANCE FROM NASA J ' ,EI 
DISC RADIUS FROM NASA 1 ' , ERR 


' . RHS 
, RZ 
' , AB 
' , RDB 
' . ABC 

,af 

= ' , BET 


, TE 
- TJ> 
, TF 






xr- . . . . 

XRN r - XD 
PRINT*. 'fcft MD t:. ' * HUH 
PRINT* i ' XE = ' ,XE 
PRINT*:* ' CALCULATED 

p e- w* r = ( j t i i ■ ? / v n< > i \ * + /» 


, PRAMDTI. 

/ Y P f.l > * 1 OO 4 A 


WAVE!. FNHH 3 S ' y X!:N 



( 


( 





c 

€ 

C 


PROGRAM' PROGIO 

s§m 



PSI ,R,CFI 


5 ~''- PRINTS. 'ENTER NOZZLE TOTAL PRESSURE IN PSIAT 7 ' > 
REAH*.PTN ' 

IF(PTN.EQ.O.O) THEN 

GOTO 400 

ENDIF 

PRINT*, 'ENTER NOZZLE TOTAL TEMPERATURE IN DEG R' ' 
READ*,TN ■ 

PRINT*, 'ENTER MEASURED THRUST IN LRF'~ 

READ*.FZ 

EF = 0,90 
10 ML = 0.1 

100 FML = (ML**4 . 0 ) * ( EF-1 » 0) + (ML**2 . 0 ) * (EF*7 . 5--20 .0 > +12 ♦ 5*EF 

MU = ML + 0.05 

FMU = <MU**4.0)*<EF-1 . 0> + < MU**2 .0>* ( EF*7 , 5-20 . 0) PI 2 . 5*Er 

FMV = FML*FMU 

IF ( FMV.GT .0.0) THEN 

ML = MU .■ ' 

GOTO 100 ' " 

ENDIF • 

. 200 MA = (MU4MD/2.0 

T' ;.Vr •••: FMA <HA*#4 .0)*<EF~1 . 0) +<MA**2,0)*<EF*7. 5-20_» 01+12 ♦ 5*EF 

FMV = FMA*FML - ' - 

. : IF ( ARS < FMA ) . LT . 0 ♦ 1 E-4 ) THEN ' 

‘MT = MA ' 

. - ELSEIF<FHV.LT. 0.0) THEN 

Mil = MA 
GOTO 200 

'ELSEIF(FMV.GT.O.O) THEN 

. . ML = MA . 

FML = <ML**4.0'i*(EF-l .0> + (ML**2.0)#(EF*7.5-20,0>U2«5*EF 
• • GOTO 200 

ENDIF >.■ 

AN = 0.162597 

P = ( 0 . 2*MT**2 . 0 ) / ( EF* ( 1 . 04 0 » 2*MT **2 . 0 ) ) 

PT = PTN* < ( 1 , 0-P ) **3 . 5 ) 

TT = TN/< 1 . 0+0 «2*MT**2 . 0) 

DT = (PT/TT ) *2.7001 18 
PI = 0.52S28*PTN 
TI = 0.8333333*TN 
DI = (PI/TT)*2. 700118 
VT = MT*49.032*(TT**0.S> 

VI = 49 . 032* ( TI **0 . 5 ) 

AT = FZ/ < 1 . 4*PI+PI-14 .7) 

PTR = PTN/14.7 
CDA= <DT*VT)/(DI*VI> 

MOOT = 0. 532063 1*PTN* AT/ (TN**0 . 5) 

MDOTA = DT*VT*AN/144 . 0 
MDOTP = ( FZ- (PI-14.7) *AT *32.2) /VI 
MDOTI = 0.5320631*PTN*AN/(TN**0 .5) 

MD = DI*VI " 

MDA = DT*UT 


aV 



w 


300 


CD = MOOT /MOOT I 

FAC FT*MT*AN*( 0 , 9193? > / ( TT**0 .5) 

r IF <HMTfr^tr^itO<IT.)'-'.THeN - V'-vC- ! - 

EF = FF r 0*001“^ • ' v 

GOTO 10 - ' 

ENDIF 


- MOOT 


+ (PI-14 . 7) *AT 
+ (PI-14. 7>*AT 


FI =(MDOT*V I 1/32.2 
FZC « MD0T*VI/32.2 
R = PT/PTN 
PSI = (1,0+0.2*MT**2.0>**0.5 
CF = FZ/(PTN*AN) 

CFA = (FZ*CD*0«5320631 )/( (TN**0.5)*MD0T) 

CFB = CD*0. 01019* ( (5.0*EF*(1.0-R**0.2857143)>**0.5)+ 
( CD*0 ♦ 5787* ( R - 14.7/PTN) )/<MT*R*PSI) 

ISPI = FZ/MDOTI 
ISPA *. FZ/MDOT 

CFI = 1.2678745 - 14.7/PTN r 


PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PR i NT* , 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*, 
PRINT*. 
PRINT* , 
PRINT*, 
PRINT*, 
PRINT*, 

PRINT*, 

PRINT*, 

PRINT*, 

PRINT*, 


'PTR = 
'EF = ' 
'MT = ' 
'AT = ' 
'AN = ' 
'PTR = 


'DT 

'FT 

'TT 

'DI 

'PI 

'TI 

'VT 

VI 


' ,PTR 
»EF 
,MT 
rAT 
, AN 
PTR 


'»4t 

' rPT 

7 »tt 

' , D I 
',PI 
' ,TI 
7 f VT 
VA - '.VI 

7 MDOT(ACT)= '.MOOT 

7 MDQTA = ' . MDOTn 

'MDOT ( I SEN > = ' , HDQTI 

7 MD = ' ,MD 

'MDA = ' , MDA 

'FFA = '.ERA 

'FZ = '.FZ 

'FI = ' ,F I 

7 R = ' ,R 

7 PSI = ',P£I 

'CD = ' ,CD 

'CD A = ' , CDA 

'CFI = ' , CFI 

'CF = '.CF 

'CFA = ' , CFA 
'CFB = ' , CFB 
'ISPI - ',ISP.I 
'ISPA = ' , I SPA 


joo G 0 T 0 5 

400 END 

END OF FILE REACHED - PRESS RCTURN 




Appendix H 


NASA Graphical Data From Ref (3) and (13) 
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Figure 13. — Effect of the ratio of specific heats of the jet upon the initial inclination of the 

jet boundary. A/,— 1.00. 
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Appendix I 


Nozzle Geometry 





15 deg Nozzle 





Appendix J 


Thrust Stand Test Data 























































































































































































































System identification data 


/ . - Chan. data *+* 

X : e 

A 7*F‘r essure 1 
A 8*Fressure 2 
A15*F'r ess 1 

A 16 *Rress 2 

A17*R1 

A13*F<2 

A19*R3 

A20*R4 

A2 1 -*R5 

A22*R6 


***•■*■* 2-24-91 

Svstem identification data 


Value Unit 


29.2 

pcig 


psig 

30 - 0 

psig 

38. 3 

psiq 

19.48 

lbf 

-17.49 

lbf 

30.54 

lbf 

5. 15 

lbf 

-1 . 13 

lbf 

-1.94 

lbf 


*** Chan. data. *** 


# Name 

Value Unit f- 

A 7* Pressure 1 

34.3 psig 

A S*Pressure 2 

43.7 psig 

-15*Pressl 

34.7 psig 

/ *Pr255 2 

44.0 psig 

l a 

19.71 lb-f 

Alfe*fc2 

-17.74 1 b-f 

A19*R3 

31-31 1 b-f 

A20*R4 

5.34 lbf 

A2 1 *R5 

-1.48 lbf 

A22*R6 

-2.43 lbf 

***** 2-24-91 

System identification 

data 

*** Chan. data **+ 

# Name 

Value Unit t 

A 7*Pressure 1 

40.6 psig 

A 8*F'ressure 2 

51.2 psig 

A15*Pressl 

40.9 psig 

A16*Press2 

51.6 psig 

A17*R1 

20.07 lbf 

A1B*R2 

-17.94 lbf 

A19*F:3 

32.41 lbf 

A20*R4 

5.64 lbf 

^ *R5 

-1.95 lbf 

■ R6 

-3.03 lbf 

***** 2—24—91 

System i dent i f i cati on 

data 


*** Chan. data *** 
# Name 


Value Uni 


A 7*Fr essure 1 

A 3 *Pr essure 2 

A15*F'ress 1 

A16*Press2 

A17*R1 

A 1 0*R2 

A19*R3 

A20*R4 

A21 *R5 

A22*R6 


**•*■** 2—24—91 

Svstem identification data 


34 . 0 

PS 

43.2 

P9 

34. 0 


43. 1 

ps 

19.73 

1 b 

-17.68 

lb 

31. 18 

lb 

5 . 33 

lb 

- 1 . 23 

lb 

-2. 15 

lb 


*** Chan. data *** 

t Name 

A 7 * Pressure 1 

A 8*F-r essure 2 

A 15 *Pressl 

A 16 *Press 2 

A 1 7 *F( 1 

A18*R2 

A19*R3 

A20*R4 

A2i *R5 

A22*R6 

*■»**•*- 2—24—91 
System identi+icat 

« *■* Chan .data *** 
# Name 


Va ! uC' 


*11 . 


• In i 1 

> 

51 ,S p " i 
4.1 . O p.si 

r , n •; i 


20 . 1 jl 
■17.98 
72. 39 
5.63 
- 1 . 72 
-2 . 80 


1 D1 

lbf 
t bf 
lbf 

x h; 

;bf 


A 7*Pr essure 1 
A 8*Pr essure 
A 15 *Pr ess 1 
A 1 6*Press2 
A 1 7 * R t 
A18*R2 
A19*R3 
A70*R4 
A21 *R5 
A22*R6 


r j a t a 

--- . 

Value Ur.it 

47. 2 

psi 

r ■ : » 

psi- 

4 2 - 0 

psi- 

52- 8 

psi < 

20, 15 

lbf 

- 1 . V8 

\ b t 

— - c„. l: 

_ \ ^ m w 4 w l 

J bi 

5.66 

lbf 

-1.99 

lbf 

-3- 10 

lbf 


** *** 2-24-91 

System identification data 


